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Abstract 

Wind  can  and  often  does  significantly  affect  impact-point  prediction  (IPP)  performance  for  thrusting/  ballistic  endoatmospheric 
projectiles.  Wind  exacerbates  the  estimation  ambiguity  between  drag  and  thrust  in  the  dynamic  model  and  induces  additional 
uncertainty  in  the  IPP  procedure.  A  tracker  accounting  for  the  wind  effect  is  presented  and  simulation  study  shows  that  it  can  be 
fully  compensated  if  the  wind  information  is  available.  An  N-point  adaptive  initialization  based  on  a  goodness-of-fit  test  and  a 
statistical  significance  test  is  introduced.  Based  on  the  multiple  interacting  multiple  model  (MIMM)  approach  developed  recently, 
the  IPP  performance  is  investigated  with  respect  to  the  total  observation  time  and  the  sensor  accuracy  in  various  wind  scenarios. 

In  each  Monte  Carlo  (MC)  run  of  the  simulation  study,  under  the  same  sensor  accuracy  and  the  same  observation  time,  the 
same  set  of  random  numbers  has  been  used  (but  different  in  different  MC  runs)  for  the  same  caliber  projectile  in  various  wind 
scenarios  to  examine  how  much  the  wind  affects  the  IPP  performance  with/without  the  exact  knowledge  of  the  wind  information. 
The  final  conclusion  is  that  with  the  wind  effect  accounted  for,  the  IPP  performance  in  the  presence  of  wind  is  practically  the 
same  as  in  its  absence. 


Index  Terms 

multiple  interacting  multiple  model  estimators,  extended  Kalman  filter,  impact  point  prediction,  wind  effect 
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2 


1.  Introduction 

Impact  point  prediction  (IPP)  of  a  thrusting/ballistic  projectile  is  significantly  affected  by  wind.  Tbe  wind  effect,  in  addition 
to  tbe  short  available  observation  time  and  the  limited  sensor  accuracy,  is  an  issue  of  high  importance  and  concern. 

The  wind  effect  causes  the  moving  projectile  to  turn  into  the  “apparent”  wind,  i.e.,  wind  causes  the  nose  of  the  moving 
projectile  to  be  pointed  into  the  wind  as  observed  from  the  projectile’s  (moving)  frame  of  reference  [6]  [18].  A  tail-  or  head¬ 
wind,  if  unaccounted  for,  would  cause  the  IPP  to  give  under-  or  over-prediction  for  the  projectile’s  range  (hence  the  term  range 
wind  effect);  a  cross  wind  would  lead  to  left-  or  right-deviation.  The  vertical  wind  is  ignored  in  this  study,  as  it  is  typically 
not  present.  Consequently,  the  projectile’s  dynamic  model  is  modified  from  [14]  [15]  to  incorporate  the  wind  effect.  Based  on 
the  modified  model,  the  sensitivity  of  the  IPP  performance  with  respect  to  the  wind  effect  is  studied. 

The  available  observation  time  is  a  key  factor  that  significantly  affects  the  IPP  performance  and  becomes  even  more  critical  in 
the  presence  of  wind  [6]  [13].  For  the  IPP  problem,  using  a  state  model  with  a  drag  coefficient  and  thrust  estimated  as  separated 
components,  there  is  an  important  observation-time  tradeoff  between  accuracy  and  implementation  of  countermeasures.  Before 
the  prediction  procedure  starts,  the  collected  measurements  are  used  to  initialize  the  estimator,  attack  the  estimation  ambiguity 
between  the  drag  coefficient  and  the  thrust  components,  and  detect  the  probability  of  the  correct  mode  converge  close  to  unity 
[15].  The  wind  effect  could  impact  each  of  these  and  evidently  deteriorate  the  overall  IPP  performance. 

Due  to  the  nonlinear  nature  of  the  dynamics  of  the  thrusting/ballistic  projectiles,  linearization  and  discretization  of  the 
nonlinear  system  is  an  approximation  of  reality.  For  a  fixed  sampling  rate,  the  approximation  quality  is  determined  by  the 
sensor  accuracy.  The  observation  uncertainty  caused  by  the  sensor  errors  (in  range  and  azimuth/elevation  angles)  affects  the 
estimation  result,  as  well  as  the  prediction  procedure  that  follows,  and  hence  the  final  IPP  performance.  Wind  uncertainty 
exacerbates  the  observation  uncertainty. 

The  present  paper  aims  to  offer  an  extensive  analysis  of  the  wind  effect  on  the  recently  developed  IPP  system  using  a 
multiple  interacting  multiple  model  (MIMM)  estimator  (with  different  initial  drag  coefficient  estimates  and  using  unbiased 
mixing^)  [15].  As  an  important  practical  issue,  in  many  practical  situations  the  wind  information  is  not  available.  To  examine 
how  much  effect  this  has  on  the  IPP  performance,  three  situations  are  considered,  i.e.,  given  perfect  knowledge  of  the  wind 
information  (denoted  concisely  as  “g.w.i.”  for  later  use),  given  no  knowledge  of  the  wind  information  (“n.w.i.”),  given  the 
wind  information  with  certain  deterministic  error  (“e.w.i.”).  In  order  to  provide  a  comprehensive  insight,  the  IPP  performance 
is  investigated  with  respect  to  the  total  observation  time  and  the  sensor  accuracy  under  various  wind  scenarios.  It  is  worth 

*The  unbiased  mixing  of  IMM  estimator  is  the  key  for  correct  estimation  of  the  extra  components  when  the  mode-matched  filters  of  IMM  estimator  are  of 
unequal  dimensions,  e.g.,  the  thrust  is  the  extra  component  in  the  present  discussion. 
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noting  that,  for  simplicity,  we  assume  a  constant  wind  during  the  estimation  and  prediction  procedure  (i.e.,  the  whole  trajectory 
period).  Altitude/location-dependent  wind  can  be  easily  dealt  with  in  the  same  manner. 

An  7V-point  adaptive  initialization  based  on  a  goodness-of-fit  test  and  a  statistical  significance  test  is  introduced.  For  each 
IMM  estimator  (associated  with  a  selected  drag  coefficient)  from  the  MIMM  estimator,  the  initialization  of  the  thrust  component 
is  of  special  importance  when  only  a  very  short  total  observation  time  is  available.  A  good  initialization  will  alleviate  the 
estimation  ambiguity  between  the  drag  coefficient  and  the  thrust.  The  A^-point  adaptive  initialization  particularly  provides  an 
early  and  reasonably  accurate  estimate  between  the  drag  coefficient  and  thrust  (as  well  as  the  kinematic  components)  and  leads 
to  a  quick  identification  of  the  correct  mode  of  the  IMM  estimator. 

The  paper  is  organized  in  following  manner.  In  Section  II,  the  wind  effect  is  presented  and  under  a  flat-Earth  assumption 
(suitable  for  short  range  projectiles),  the  dynamic  model  in  the  presence  of  wind  is  modified  from  [15].  The  corresponding 
discretized  form  of  the  modified  dynamic  model  and  the  discrete-time  measurement  equations  are  also  presented  here.  The 
MIMM  estimator  and  the  IPP  procedure  are  described  in  Section  III.  The  parameter  setting,  design  of  the  MIMM  estimator 
and  the  7V-point  adaptive  initialization  approach  are  presented  in  Section  IV.  Then  in  the  presence  of  different  wind  (strength 
and  direction)  conditions,  the  IPP  performance  is  investigated  in  various  total  observation  time  and  sensor  accuracy  scenarios 
by  simulation  in  Section  V.  Conclusions  are  presented  in  Section  VI. 

II.  Dynamic  Model  and  Measurement  Model  in  the  Presence  oe  Wind 

The  wind  effect,  which  worsens  the  estimation  ambiguity  between  the  drag  coefficient  and  thrust  [15]  and  thus  presents 
an  additional  challenge  to  the  IPP,  must  be  carefully  accounted  for  [3].  The  sensitivity  of  the  IPP  performance  to  the  drag 
coefficient  estimate,  which  is  significantly  affected  by  the  relative  velocity  of  the  projectile  (even  one  with  known  shape)  with 
respect  to  air,  necessitates  the  quantification  of  the  wind  effect  for  the  purpose  of  IPP. 

The  contribution  of  the  wind  comprises  the  range  wind,  cross  wind  and  vertical  wind  effects.  We  ignore  the  last,  since  it  is 
generally  small;  but  the  techniques  we  present  could  be  augmented  to  account  for  it.  A  range  (head/tail)  wind  will  push  back 
or  forward  the  moving  object  and  a  cross  wind  causes  the  moving  object  deviate  to  the  side.  These  wind  components  together 
are  considered  as  the  true  wind  velocity  and  are  combined  with  the  target  velocity  with  respect  to  the  ground  to  yield  the 
so-called  “apparent  wind”.  The  moving  target  turns  into  the  apparent  wind,  i.e.,  its  nose  is  pointed  into  the  wind  (or  “upwind”) 
while  its  tail  pointing  “downwind”  [6]  [18].  With  the  wind  effect  shown  in  Fig.  1,  we  have 

Va=V^,-V  (1) 

where  v  is  the  target  velocity  with  respect  to  the  ground,  v^,  is  the  wind  velocity  and  is  the  apparent  wind  velocity. 
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(Projectile) 


Fig.  1:  The  wind  effect  (v  is  the  target’s  ground  velocity,  v^,  is  the  wind  velocity  and  is  the  apparent  wind  velocity). 


Thus,  the  drag  and  the  thrust  in  the  following  dynamic  model  should  be  aligned  with  the  direction  of  the  projectile,  which  is 
aligned  with  the  apparent  wind.  In  the  absence  of  information  about  how  fast  the  projectiles  align  themselves  with  the  apparent 
wind  following  launch,  it  is  assumed  that  this  happens  instantaneously.  During  thrusting  phase  this  is  especially  important,  as 
due  to  wind  thrust  may  thus  not  be  aligned  with  the  direction  of  travel. 

For  simplicity  of  analysis,  we  assume  a  constant  wind  during  the  estimation  and  prediction  procedure  (i.e.,  the  whole 
trajectory  period).  The  wind  environment  where  the  instantaneous  projectile  alignment  hypothesis  is  reasonable  and  acceptable 
(e.g.,  altitude/location-dependent  wind  with  slow  variation  or  time-variant  wind  that  is  piecewise  constant  with  respect  to 
altitude)  can  be  dealt  with  in  the  same  manner. 

In  the  presence/absence  of  wind,  the  trajectory  of  a  thrusting  projectile,  from  launch  to  impact,  can  be  divided  into  two 
phases;  thrusting  and  ballistic.  It  is  a  natural  choice  that  we  use  an  IMM  estimator  with  a  thrust  mode  (TM)  and  a  ballistic 
mode  (BM)  to  match  these  phases.  The  state  vector  for  the  corresponding  TM  is 

x{t)  =  [x{t)  y(t)  z{t)  x{t)  y{t)  z{t)  a{t)  r(t)]'  (2) 


where  a{t)  is  the  drag  coefficient  and  r(f)  is  the  thrust.  The  state  vector  of  the  BM  is  the  same  as  above  but  excludes  the 
thrust  component.  The  time  arguments  will  be  omitted  where  no  ambiguity  is  caused. 

Under  Hat-Earth  assumption,  the  dynamic  model  in  the  presence  of  wind  modified  from  [15]  (the  corresponding  modifications 
are  indicated  by  the  sub-script  “a”,  which  implies  the  use  of  “apparent”  wind  information)  can  be  written  as  follows 


~x' 

'Xa' 

y 

II 

Va 

Oirn 

Va 

-Z. 

-Za  - 

-Za- 

S  +  vi 


(3) 


and 


a  =  V2 

(4) 

T  =  Ds 

(5) 
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where 

•  The  target  velocities  with  respect  to  air  are  =  x  —  x^,  i/a  =  y  —  Vw  Za  =  z  —  Zyj. 

•  The  first  term  on  the  right  side  of  (3)  is  the  specific  thrust  in  the  x,  y,  and  z  directions  with  accounting  for  the  wind 

effect.  For  the  ballistic  phase,  the  thrust  is  zero. 

•  The  second  term  is  the  drag  part,  which  is  related  to  the  altitude  and  target  velocity  with  respect  to  air. 

•  Va  is  the  magnitude  of  the  velocity  [  Xa  Va  Za  ]',  i.e.,  the  speed  (in  m/s). 

•  a  is  the  drag  coefficient  (in  m^/kg)  at  subsonic  speed  and  r  is  the  thrust  (in  m/s^).  It  is  known  that  the  drag  coefficient 

varies  significantly  with  the  Mach  number  regime:  subsonic,  transonic  and  supersonic.  This  will  be  accounted  for  by  the 

am  as  below. 

•  am  is  the  (dimensionless)  Mach  number-dependent  drag  coefficient  multiplier,  which  is  approximated  by  a  cubic  spline 
curve  shown  in  Fig.  2.  The  cubic  spline  curve  is  obtained  by  an  interpolation  process  shown  in  detail  in  Appendix  A. 

•  Da  =  —  ^  yyjjgj-g  is  the  air  density  (in  kg/m^)  at  altitude  z  (in  m)  and  c  is  the  air  density  constant 

(in  m-i)  [12]. 

•  g=[0  0  — pj'is  the  standard  gravity  vector,  g  is  the  standard  acceleration  due  to  gravity  at  sea  level,  assumed  to  be 
the  same  throughout  the  trajectory,  with  value  9.812  m/s^. 

•  vi,  V2,  and  F3  are  assumed  to  be  continuous-time  zero-mean  white  Gaussian  noises.  The  drag  coefficient  and  thrust 
acceleration  are  thus  modeled  as  Wiener  processes  with  slow  variation  [1]. 

Combining  the  dynamic  equations  (3)-(5),  we  have  the  following  compact  form  (accounting  for  the  wind  effect) 

x(f)  =  /u,[x(f),  x„(f)]  -b  v(t)  (6) 


where 


/„[x(f),  x^(f)] 


with  the  wind  effect  compensation  vector 


x{t) 

y{t) 

z{t) 

+  a{t)am(t)Da{t)Xait) 

+  a{t)am{t)Da{t)ya{t) 

+  a{t)am{t)Da{t)za{t)  -  g 
0 
0 


(7) 


Xu.(f)  =  [0  0  0  Xa,{t)  yw{t)  Zyj{t)  0  Oj' 


(8) 
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Fig.  2;  The  cubic  spline  approximation  of  the  Mach  number-dependent  drag  coefficient  multiplier  (for  a  “sharp  nose” 

projectile  [17]). 


and  the  Wiener  process  vector 

=  V2{t)  F3(t)]'  (9) 

Note  that  /«(•)  and  /r(-)  are  both  equal  to  zero,  as  we  assume  the  drag  coefficient  and  thrust  are  Wiener  processes. 

The  state  vector  equation  is  discretized  by  a  second  order  Taylor  expansion^  [7].  Using  the  discrete  time  notation  x(fc)  = 
x(f)|t=feT  and  x(fc  +  1)  =  x(f  +  T)\t=kT  (the  same  for  we  have  the  following  discretized  continuous-time  expression 

[15] 

x(fc  -f  1)  =  K{k)  +  f^[x{k),Xyj{k)]T  +  A^{k)fyj[x{k),K^{k)]'^  +  v{k)  (10) 

where  Ayj{k)  is  the  Jacobian  of  (7)  evaluated  at  x(fc)  (with  known  wind  information  x^{k))  and  i2{k)  is  the  discretized 
continuous-time  (white  Gaussian)  process  noise  for  the  sampling  interval  T.  Based  on  the  assumption  that  a  is  nearly  constant 
and  Da  is  related  to  both  z  and  Va,  the  detailed  form  of  A^j  here  is  obtained  by  replacing  in  the  Jacobian  in  [15]  the  target 
velocities  (with  respect  to  ground)  with  the  corresponding  target  velocities  with  respect  to  air. 

^It  has  been  reported  that  for  a  short  range  ballistic  projectile  the  first  order  Taylor  expansion  is  sufficient.  However,  for  the  thrusting  projectile  with 
unknown  drag  coefficient  and  unknown  thrust  -  as  in  the  present  study  -  the  second  order  Taylor  expansion  is  used  to  compensate  for  the  nonlinearity  made 
obscure  by  the  unknowns.  The  methodology  can  be  used  for  long  range  scenario,  in  which  case  a  second  order  Taylor  expansion  is  necessary  to  avoid  bias 


[2]. 
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The  corresponding  covariance  matrix  of  the  discretized  process  noise  is 


Q  = 


r  7-  T 

-^13 


Th 

0lx6 


Qv  Ogxi  Oexi 

Tqa  0 


(11) 


0lx6  0  Tq-r 

where  I3  is  the  3x3  identity  matrix  and  the  continuous  time  process  noise  “intensities”  q^,  qa  and  q^  are  the  corresponding 
power  spectral  densities  (PSD). 

Assuming  the  sensor  is  located  at  (x*  ys  Zg  ),  the  measurements  in  spherical  coordinates  are 


=  r  +  Wr  =  a/(x  -  XaY  +  {v  -  VsY  +  {z  -  ZgY  +  Wr 
-1  f  y-Vs 


0m  =  0  +  wg  =  tan 


X  — 


■  wg 


Em  =  e  +  it'e  =  tan 


-1 


2  —  2, 


(12) 

(13) 

(14) 


\/{x-  XsY  +  {y-  VsY  +  {z-  ZaY  j 

where  r,  9  and  e  are  the  independent  true  range,  azimuth  and  elevation  components,  respectively.  Wr,  wg  and  denote  the 
corresponding  zero-mean  white  Gaussian  measurement  noises  with  standard  deviations  (SD)  ar,  erg  and  respectively. 

An  unbiased  measurement  conversion  from  spherical  to  Cartesian  coordinates  was  presented  [2],  so  that  the  IPP  problem 
can  be  described  entirely  in  Cartesian  coordinates: 


z(fc)  =  Hx{k)  +  w(fc) 


(15) 


where  H  =  [I^  0],  w(fc)  is  the  equivalent  measurement  noise  vector  in  Cartesian  coordinates  obtained  from  the  unbiased 
measurement  conversion  and  R{k)  is  the  corresponding  equivalent  state-dependent  covariance  matrix. 


III.  Multiple-IMM  Estimator  for  IPP 

During  the  thrust  mode  (TM),  the  drag  parameter  (drag  coefficient)  and  thrust  parameter  are  separate  state  components  as 
shown  in  (2).  However,  the  drag  force  and  thrust  force  are  acting  simultaneously  (see  equation  (3))  and  the  IMM  estimator 
has  difficulty  distinguishing  between  them  if  the  initial  uncertainty  in  the  drag  coefficient  is  large.  A  sudden  decrease  in  the 
drag  coefficient  estimate  may  trigger  an  increase  in  the  thrust  estimate.  In  the  presence  of  wind,  the  wind-induced  motion 
makes  correctly  unsnarling  these  important  components  (drag  coefficient  and  thrust)  increasingly  obdurate.  The  sensitivity  of 
the  estimation  to  the  initial  drag  coefficient  estimate  necessitates  the  use  of  an  MIMM  estimator  to  overcome  this  “marginal 
observability”  problem  [1]. 

The  procedure  starts  by  establishing  a  set  of  the  L  IMM  estimators,  each  with  an  appropriate  set  of  modes  (TM  and  BM 
in  present  case)  to  describe  the  system  behavior.  Each  IMM  estimator  will  be  initialized  with  a  different  value  of  the  drag 
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coefficient  estimate  with  a  suitable  initial  SD.  The  initial  SD  is  taken  equal  to  25%  of  the  initial  estimate  of  the  drag  coefficient. 
The  filtering  parameter  setting  is  discussed  in  Section  IV. 

In  order  to  select  the  best  initial  drag  coefficient  estimate,  we  need  to  determine  the  most  likely  among  L  IMM  estimators, 
during  the  observation  period.  The  likelihood  function  of  IMM  estimator  I  for  the  time  interval  [/cq,  K]  is 

K 

^  n  A,  (fc)  1  =  1,...,  L  (16) 

k—ko 

where  A;  (k)  is  the  likelihood  function  of  Ah  IMM  estimator,  with  r  modes,  at  time  k,  given  by  [2] 

M{k)=  (17) 

where  A;(fc)  is  the  likelihood  function  of  mode  i  of  IMM  estimator  I  at  time  k  and  /r;(fc|fc—  1)  is  the  predicted  mode  probability 
for  mode  i  of  IMM  estimator  1.  The  mode  likelihood  is  [1] 


A]{k)  =  p  [z{k)\Mi{k),  xO*(fc  -  l|fc  -  l),P°\k  -  l\k  -  1)] 

=  Af[2{ky,zl{k\k-l),Stik)] 


(18) 


where  A/"  [  •  ]  is  the  Gaussian  probability  density  function,  zj  is  the  predicted  measurement  and  is  the  innovation  covariance 
in  mode  i  of  IMM  estimator  1.  The  predicted  mode  probability  can  be  written  as 

/rKA:|fc-l)= 


i=i 


i=i 


(19) 


where  pj  [k  —  1)  is  the  mode  probability  of  mode  j  of  IMM  estimator  I  at  time  k—1  and  pji  is  the  transition  probability  from 
mode  j  to  mode  i  over  one  time  interval. 

Using  (18)  and  (19),  one  can  rewrite  (17)  as 


Then 


i=i  i=i 


I*  =  argmax  =  argmin  —  InAj'”’^ 


(20) 


(21) 


selects  the  best  IMM  estimator.  This  yields  the  best  initial  estimate  of  the  drag  coefficient. 

Based  on  the  MIMM  estimator,  we  can  choose  at  the  end  of  the  observation  period  the  most  likely  initial  drag  coefficient 
of  the  most  probable  mode  from  the  best  IMM  estimator  (21).^  Then  a  numerical  open-loop  nonlinear  predictor  (we  use  the 

^In  practice,  the  thrusting  period  is  relatively  short  and  by  the  end  of  the  observation  period  the  BM  should  be  the  dominant  one  in  the  IMM.  Otherwise, 


unless  one  were  clairvoyantly  to  know  the  burnout  time,  one  cannot  make  a  meaningful  IPP. 
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4th  order  Runge-Kutta  (RK)  method^  [11],  incorporating  the  wind  effect  if  the  wind  information  is  available)  is  employed  to 
predict  the  trajectory  down  to  its  impact  point.  At  the  same  time,  the  corresponding  covariance  is  also  predicted  to  the  impact 
point  using  a  zero-gain  (open  loop)  EKF  covariance  equation.  Then  the  99%  chi-square  probability  region  ellipse  for  the  true 
impact  point  is  [1] 


-1 


[  ^ip  ^pd  Uip  Vpd  ]  Ppfi 


^ip  ^pd 

^ip  Vpd 


=  Xl  (99%) 


(22) 


where  {xip  yip)  is  the  true  impact  point,  {xpd  Vpd)  is  the  predicted  impact  point,  Ppd  is  the  corresponding  predicted  covariance 
matrix,  and  xi  (99%)  denotes  the  99%  point  on  the  chi-square  cumulative  distribution  function  with  two  degrees  of  freedom 
[1].  This  can  be  used  as  the  “warning  zone”  in  practical  situations. 

In  practical  situations  the  ellipse  is  centered  at  the  predicted  impact  point  and  one  can  check  whether  the  true  impact  point 
(if  available)  is  within  the  99%  probability  ellipse  [15].  In  Monte  Carlo  simulations,  as  in  present  study,  we  evaluate  the 
(equivalent)  converse;  whether  the  predicted  impact  point  falls  into  the  ellipse  centered  at  the  true  impact  point  [12]. 

The  IPP  procedure  based  on  the  MIMM  estimator  is  shown  in  Fig.  3. 


Fig.  3;  The  IPP  procedure  based  on  the  MIMM  estimator. 


^Using  the  second-order  open-loop  EKF  directly  with  a  small  iteration  time  interval  (0.1  s)  for  the  prediction  shows  IPP  performance  practically  the  same 


as  the  RK  method. 
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IV.  MIMM  Design  Parameter  Selection 

The  MIMM  estimator  using  unbiased  mixing  with  different  initial  drag  coefficient  estimates  is  used  [15].  The  unbiased 
mixing  is  necessary  because  of  the  unequal  dimension  state  vectors  in  the  two  modes  (TM  and  BM). 

A.  MIMM  for  selection  of  the  best  initial  drag  coefficient  estimates 

The  MIMM  estimator  is  chosen  to  consist  of  T  =  4  IMM  estimators  with  initial  drag  coefficient  estimates  a*(0),  i  = 
1,  2,  3, 4  as  0.18,  0.13,  0.065  and  0.03  (m^/kg),  respectively.^  The  initial  SD  of  each  drag  coefficient  estimate  is  25%  of  the 
corresponding  initial  estimate  of  the  drag  coefficient. 

If  the  projectile  library  with  the  truth  of  the  corresponding  drag  coefficient  information  is  provided,  the  selection  of  the  best 
initial  estimate  of  the  drag  coefficient  actually  indicates  the  projectile  identification  [4]  [12].  However,  the  wind  effect  could 
blur  the  identification. 


B.  N-point  adaptive  initialization 

The  initialization  is  crucial  to  the  accuracy  of  the  state  estimates.  Due  to  the  sensitivity  of  the  estimation  on  the  drag 
coefficient,  which  is  velocity  (Mach  numberj-dependent  and  “marginally  distinguishable”  from  the  thrust  estimate,  an  iV-point 
adaptive  initialization  method  is  used  here  (the  number  N  is  discussed  later).  This  method  is  based  on  the  polynomial  fitting  of  a 
set  of  noisy  position  measurements  [1],  with  the  polynomial  order  adaptively  chosen  based  on  the  corresponding  goodness-of-fit 
error.  Particularly,  it  can  give  a  good  guideline  for  how  to  initialize  the  key  thrust  component. 

The  iV-point  fitting  of  a  polynomial  of  order  n  is  done  as  follows.  The  position  measurements  for  the  three  Cartesian 
coordinates  are 


where  k  =  1,2, ...  ,N  and 


z(fc) 


=  Ha{k)a  +  w{k) 


1 

o 

o 

1 _ 

A 

o 

o 

Ry 

1 

o 

o 

_ 1 

.Hz. 

+  w{k) 


h{k) 


1  ffc  ••• 


jtk) 

n\ 


(23) 


(24) 


with  tk  is  the  sampling  time  and  the  parameter  vectors 


=  [  ttio  an  •  •  •  Qin  ]  i  =  X,  y,  z  (25) 

^This  grid  of  values  was  chosen  based  on  the  tracking  results  in  the  absence  of  wind.  It  seems  that  the  60  mm  projectile  has  a  somewhat  different  drag 


coefficient  in  the  presence  of  wind  (it  turns  somewhat  slowly  into  the  wind,  perhaps  the  result  of  a  larger  moment  of  inertia). 
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contain  the  coefficients  of  the  polynomials  (one  for  each  coordinate). 

For  the  time  covered  by  N  points,  we  get  the  estimate  for  the  parameter  vector  a  (of  dimension  Ua  =  3(n  +  1))  as 


H^'  (R^)  ^  *  Hf  (R^) 


-1 


jN'  {uN\-^ 


with  the  corresponding  covariance  matrix 


where 


,{l) 


z^  = 


z(7V) 


is  the  stacked  vector  of  measurements  (of  dimension  iN  x  1), 


= 


H.{1) 


is  the  stacked  measurement  matrix  (of  dimension  3N  x  3(n  +  1)),  and 

-R{1)  •••  0 

R^  =  :  : 

.0  •••  R{N)} 

is  the  block  diagonal  covariance  matrix  of  measurement  noise  (of  dimension  3N  x  3iV)  for  the  fitting  interval. 
The  goodness-of-fit  error  has  the  following  chi-square  distribution  [1] 

A 


Jn  =  -  iTf  a]'  (rY  '  [z^  -  ~  xIn-. 


(26) 


(27) 


(28) 


(29) 


(30) 


(31) 


Using  as  acceptance  region  for  (31)  its  95%  probability  region  (one-sided)  and  a  component  statistical  significance  test  (a 
Gaussian  test  for  the  absolute  value  with  95%  probability  region,  i.e.,  two-sided),  we  can  adaptively  choose  the  best  order 
n*,  which  avoids  both  “over-parameterizations”  and  “under-parameterizations”  [1].  This  is  done  by  starting  with  n  =  1  and 
increasing  it  until:  (i)  the  test  statistic  (31)  falls  below  the  95%  probability  threshold  and  (ii)  the  magnitude  of  at  least  one  of 
the  components  (x,  y  and  z)  is  statistically  significant  with  the  threshold  Q{97.5%),  which  amounts  to  cutting  the  upper  and 
lower  tails  of  2.5%. 

The  7V-point  adaptive  initialization  used  =  12  (about  1  s  of  data).  As  an  illustrative  example.  Table  I  shows  the  results 
of  the  polynomial  fitting  of  an  observation  sequence  of  the  trajectory  60C7H7®  with  constant  velocity  (CV,  n  =  1),  constant 
acceleration  (CA,  n  =  2)  and  constant  jerk  (CJ,  n  =  3)  models  in  one  of  the  MC  runs.  For  this  particular  case,  the  CA 

®This  stands  for  caliber  60  mm  projectile  in  the  presence  of  7  m/s  crosswind  and  7  m/s  headwind. 
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(n*  =  2)  model  is  chosen;  the  fitting  error  is  =  44.0  <  X3iv-na(^®^)  “  z-component  estimate  magnitude  is 

statistically  significant:  2.3  >  Q(97.5%)  =  1.96.  The  test  statistics  that  yielded  the  choice  n*  =  2  are  in  boldface. 

Remarks:  Initialization  of  Thrust 

In  the  present  simulation  study,  the  CV  (n  =  1)  model  and  the  CA  (n  =  2)  model  may  be  selected.  For  n*  =  2,  the 

acceleration  estimate  is  used  (after  subtracting  the  gravity  acceleration  and  drag  vector)  to  obtain  the  initial  estimate  of  the 

thrust,  t(0);  for  n*  =  1  (low  thrust  case),  then  f(0)  =  g  with  SD  g/A. 

TABLE  I:  Fitting  of  CV,  CA  and  CJ  models  for  trajectory  60C7H7  (with  cr^.  =  10  m  and  ag  =  =  5mrad) 


Model 

CV  (n  = 

1) 

CA  (n  = 

2) 

CJ  (n  =  3) 

Component 

X 

y 

z 

X 

y 

z 

X 

y 

z 

O-iO 

0.9 

-5.6 

7.8 

5.8 

-17.2 

-14.4 

10.3 

-20.9 

-11.9 

VPiO 

6.8 

8.3 

10.4 

9.3 

11.3 

14.2 

11.0 

13.3 

16.7 

\o-io\/ '/Pm 

0.1 

0.7 

0.8 

0.6 

1.5 

1.0 

0.9 

1.6 

0.7 

an 

119.9 

117.2 

155.9 

91.4 

187.2 

289.2 

29.3 

238.6 

253.1 

10.4 

12.6 

15.8 

39.1 

47.4 

59.4 

90.2 

108.9 

136.7 

|<iil|/ VFil 

11.5 

9.3 

9.8 

2.3 

4.0 

4.9 

0.3 

2.2 

1.9 

ai2 

51.2 

-127.4 

-242.2 

344.6 

-369.4 

-69.0 

VP^ 

68.3 

82.7 

103.7 

389.4 

471.8 

591.3 

\ai2\l\fpv2 

0.8 

1.5 

2.3 

0.9 

0.8 

0.1 

a-iz 

-531.9 

437.7 

-316.9 

694.7 

842.9 

1055.6 

|ai3l/\/-Pi3 

0.8 

0.5 

0.3 

Jn 

52.9 

44.0 

43.3 

xiiV  — na  (95%) 

52.0 

50.7 

49.5 

C.  Selection  of  Process  Noise  Intensities 

In  the  dynamic  equations,  we  actually  assume  a  nearly  constant  velocity  model  (continuous  time  white  noise  acceleration 
—  CWNA  [1])  for  the  kinematic  components  and  a  Wiener  process  with  slow  variation  for  the  drag  coefficient  and  thrust.  In 
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order  to  satisfy  this  model’s  assumptions,  we  choose  small  process  noise  “intensities”  (power  spectral  densities  —  PSD)  q^, 
Qa,  and  Qt  as  follows. 

The  process  noise  induced  root  mean  square  (RMS)  rate  of  change  in  velocity  over  an  interval  A  is 

dv  =  [(m/s)/s]  (32) 


For  the  drag  coefficient  the  rate  of  change  is 


For  the  thrust  the  rate  of  change  is 


da  =  [(m^/kg)/s] 


A  a/^tA  /  2\  /  1 

dr  =  ^  [(m/s  )/s] 


Then  we  have  the  following  estimator  design  procedure  for  selection  of  the  process  noise  PSD  based  on  the  corresponding 
RMS  rates  of  change 

=  dlA  (m^/s^)  (35) 

qa  =  dlA  (m‘‘/(kg^-s))  (36) 

qr  =  dlA  (mVs^)  (37) 

The  process  noise  intensities  are  chosen  based  on  the  process  noise  induced  RMS  change  (in  velocity/ drag  coefficient/ thrust) 
over  an  interval  of  A  =  1  s  as  shown  in  Table  II.  Note  that  dy  is  chosen  differently  for  the  four  initial  drag  coefficient  estimates 

d*(0),  i  =  1,2, 3, 4. 

TABLE  II:  The  RMS  change  rate  due  to  the  process  noise 

dy  [(m/s)/s] 


Filter 

1 _ 

al(0) 

&\0) 

a3(0) 

&\0) 

da  [(m^/kgVs] 

dr  [(m/s^)/s] 

IMM  (TM) 

3.5 

2.5 

2 

1 

0.12dT0)/s 

0.25f(0)/s 

IMM  (BM) 

1.6 

1.4 

1.2 

0.9 

0.10dT0)/s 

N/A 

D.  IMM  Mode  Transition  Probability  Matrix 

The  Markov  chain  probability  transition  matrix  for  the  two-mode  IMM  is 

Pit  I-Pn" 

1  —  P22  P22  - 


(38) 
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where  the  elements  of  the  matrix  are  obtained  based  on  the  mean  sojourn  time  (MST)  [1],  si  and  S2,  for  the  TM  {i  =  1)  and 
the  BM  {i  =  2)  respectively 

^  i  =  1, 2  (39) 


T 

Pii  =  i - 

Si 


and  T  =  0.1s  is  the  (fixed)  sampling  interval.  Setting  the  MST  as  si  =  2s  and  S2  =  50s  for  TM  and  BM  respectively,  we 
get  the  mode  transition  probability  matrix 

r  0.950  0.050] 

(40) 


0.950  0.050 
L  0.002  0.998 


with  initial  mode  probability  vector  [0.90  0.10] 


V.  Simulation  Results 

Three  categories  of  thrusting/ballistic  projectiles  of  different  calibers  are  considered:  60  mm,  81mm  and  120  mm.  In  all, 
21  trajectories  were  generated  using  [16]^  with  a  Hat-Earth  model  with  quadrant  elevation  45°  (i.e.,  aimpoint  is  NE).  The 
thrusting/ballistic  projectiles  are  actually  rocket-assisted  with  initial  velocity  around  250  m/s  (which  varies  across  different 
calibers).  Eor  each  category,  there  are  trajectories  with  no  wind  present,  labeled  as  “WO”;  with  5  m/s  crosswind  only,  labeled 
as  “C5”;  with  10 m/s  crosswind  only,  labeled  as  “CIO”;  with  5 m/s  crosswind  and  5 m/s  headwind,  labeled  as  “C5H5”;  with 
7  m/s  crosswind  and  7  m/s  headwind,  labeled  as  “C7H7”;  with  5  m/s  crosswind  and  5  m/s  tailwind,  labeled  as  “C5T5”;  with 
7  m/s  crosswind  and  7  m/s  tailwind,  labeled  as  “C7T7”.  Note  that  the  headwind  blows  against  the  motion  of  the  projectiles 
while  the  tailwind  blows  in  the  travel  direction  of  the  projectiles  [18].  Measurements  were  obtained  with  no  missed  detections 
and  no  false  alarms. 

The  parameters  of  the  trajectories  of  the  different  caliber  projectiles  considered,  namely,  the  projectile  range  to  impact  (i?*), 
the  impact  time  (T/^),  the  sampling  interval  (T),  the  maximum  ground  speed  the  ground  speed  at  impact  point  (V]p  and 

the  apogee  altitude  are  summarized  in  Table  III.  The  launch  point  of  each  projectile  was  at  the  origin  of  coordinates. 

The  sensor  location  was  (  5000  4000  0 )  m. 

As  shown  in  Table  IV,  four  different  sensor  accuracy  settings  are  used  and  labeled  as  Case  1  (good  sensor  accuracy).  Case 
2  (poor  angle  accuracy).  Case  3  (poor  range  accuracy)  and  Case  4  (the  best  sensor  accuracy). 

With  100  Monte  Carlo  (MC)  runs  for  each  sensor  accuracy,  we  examine  the  IPP  performance  with  given  wind  information 
(g.w.i.)  and  with  no  wind  information  provided  (n.w.i.)  for  various  total  observation  time  percentages  of  the  whole  trajectory 
(denoted  as  “OT”  for  conciseness):  15%,  20%,  25%,  30%,  40%  and  50%.  In  particular,  under  the  same  sensor  accuracy  and 

^This  is  a  high  fidelity  simulation  for  rocket-assisted  projectiles,  i.e.,  with  initial  “catapulted”  velocity.  The  trajectories  were  generated  based  on  models  we 


cannot  control  and  do  not  know. 
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TABLE  III:  The  parameters  of  the  trajectories  of  the  different  caliber  projectiles  considered 


Caliber 

(m) 

r(s) 

^max 

^apg  (^) 

60  mm 

3600  ±  too 

30  ±  1 

0.1 

240  ±3 

155  ±3 

1085  ±  5 

81  mm 

5700  ±  150 

37  ±  1 

0.1 

295  ±3 

200  ±3 

1680 ± 10 

120  mm 

7000  ±  250 

41  ±  1 

0.1 

315  ±5 

230  ±5 

2055  ±  10 

TABLE  IV:  Sensor  accuracy  settings 


(Tr 

<70 

(^e 

Case  1 

10m 

5  mrad 

5  mrad 

Case  2 

10m 

lOmrad 

10  mrad 

Case  3 

25  m 

5  mrad 

5  mrad 

Case  4 

5  m 

3  mrad 

3  mrad 

the  same  OT  setting,  in  each  MC  run,  the  same  random  number  set  (but  different  in  different  MC  runs,  of  course)  has  been 
used  for  different  trajectories  (W0/C5/C10/C5H5/C7H7/C5T5/C7T7)  to  examine  how  significantly  the  wind  affects  the  IPP 
performance. 

Three  different  situations  concerning  the  knowledge  of  the  wind  information  are  considered:  the  “g.w.i.”  (perfect  knowledge), 
the  “n.w.i”  (no  knowledge)  and  the  “e.w.i”  (an  inaccurate  knowledge  of  30%  positive  deterministic  error,  i.e.,  v^,  +  = 

Vuj  +  O.Svu,,  where  v;,  is  a  deterministic  bias).  Note  that  all  the  IPP  results  are  by  default  obtained  with  “g.w.i.”  unless  “n.w.i” 
or  “e.w.i”  is  indicated.  The  indicating  term  “g.w.i.”  is  omitted  in  Eigures,  Tables  and  discussions  if  there  is  no  ambiguity. 

A  sample  estimated  trajectory,  indicating  the  portions  corresponding  to  the  MIMM  filtering  and  to  the  IPP,  including  the 
estimated  burnout  point  (BoP)  and  the  uncertainty  ellipse  centered  at  the  true  impact  point,  is  presented  in  Pig.  4. 

A.  IPP  RMS  Error 

The  IPP  root  mean  square  (RMS)  errors  of  the  trajectories  of  caliber  60  mm,  81mm  and  120  mm  (in  various  scenarios)  are 
shown  in  Tables  V,  VI  and  VII,  respectively.  The  g.w.i.,  n.w.i.  and  e.w.i.  IPP  RMS  errors  for  various  OTs  and  various  sensor 
accuracies  are  evaluated  and  compared. 

Prom  Tables  V-VII,  it  can  be  seen  that,  generally,  as  the  wind  becomes  stronger,  the  difference  between  the  g.w.i.  IPP 
RMS  error  and  the  n.w.i.  IPP  RMS  error  will  increase,  especially  for  smaller  OTs;  as  the  OT  increases,  the  influence  of  the 


wind  on  the  IPP  performance  becomes  smaller.  As  the  sensor  accuracy  becomes  worse,  the  IPP  performance  is  gradually 
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Caliber120C7T7 


_ 

MIMM  filtering 

Estimated  BoP 

IPP  procedure 

□ 

Uncertainty  ellipse 

V 

IPP 

A 

True  impact  point 

m 

Sensor 

X-axis(m) 


Fig.  4;  A  sample  estimated  trajectory  with  IPP  uncertainty  ellipse  centered  at  true  impact  point. 


degraded.  Particularly,  a  sensor  error  that  causes  more  uncertainty  along  the  travel  direction  of  the  projectiles,  compared  with 
the  sensor  error  that  causes  more  uncertainty  perpendicular  to  travel  direction  of  the  projectiles,  has  more  influence  on  the  IPP 
performance.  In  the  presence  of  a  weak  wind  (C5),  the  IPP  RMS  error  (g.w.i.)  is  practically  the  same  as  the  IPP  RMS  error 
(e.w.i.).  The  knowledge  of  the  wind  information  with  a  small  deterministic  bias  has  small  effect  on  the  IPP  performance. 

We  can  get  an  IPP  RMS  error  of  more  than  500  m  for  the  60  mm  projectile  (Table  V;  Case  2,  OT=15%),  more  than  650  m 
for  the  81mm  (Table  VI:  Case  3,  OT=15%)  and  more  than  750  m  (even  900  m)  for  the  120  mm  (Table  VII:  Case  2  and  Case 
3,  OT=15%).  These  are  unacceptable  errors.*  A  strong  wind  effect,  combined  with  a  poor  sensor  accuracy,  can  lead  to  an 
unacceptable  IPP  RMS  error  when  only  a  very  short  total  observation  time  is  available.  On  the  other  hand,  a  good  sensor 
accuracy  in  the  presence  of  known  wind  (even  with  a  small  bias)  always  yields  a  good  IPP  performance. 

In  Tables  V-VII,  the  columns  with  the  wind  information  provided  (g.w.i.)  show  that,  given  a  sufficient  observation  time  (say, 
OT=  25%-50%)  to  correctly  select  the  best  IMM  from  the  MIMM  estimator  (i.e.,  to  overcome  the  “marginal  observability” 
problem  between  the  drag  coefficient  and  thrust),  the  IPP  RMS  errors  are  practically  the  same  as  without  wind.^  This  implies 

that  the  wind  effect  can  be  fully  compensated  if  the  wind  information  is  correctly  provided. 

*In  practical  situations,  a  minimum  acceptable  IPP  RMS  error  is  expected  to  be  within  10%  of  the  projectile  range.  The  IPP  RMS  errors  larger  than  10% 

of  the  projectile  range  are  in  italic  font  in  Tables  V-VII  . 

^This  is  highlighted  by  the  IPP  RMS  errors  shown  in  bold  font  in  Tables  V-VII  . 
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TABLE  V:  IPP  RMS  errors  (in  m,  100  MC  runs)  for  various  OTs  (shown  in  the  hrst  column)  and  various  sensor  accuracies, 

caliber  60  mm:  italic,  see  footnote  8;  bold,  see  footnote  9. 


Case  1 

60W0 

60C5 

60C10 

60C5H5 

60C7H7 

60C5T5 

60C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

417.4 

419.5 

420.7 

419.8 

441.4 

426.7 

440.1 

452.5 

461 

459.6 

458.4 

464.3 

465.3 

417.2 

426.8 

418.5 

408.8 

409.6 

410.9 

20% 

152.1 

150.8 

154.7 

151.3 

150.4 

169 

157.3 

145.6 

129.5 

138.7 

142.1 

136.4 

130.2 

158.2 

182 

180.1 

157.7 

191.5 

195.2 

25% 

127 

127.1 

132.8 

127.5 

126.9 

145.2 

125.8 

122.6 

111.8 

97.4 

120.7 

119 

91.2 

127.7 

142.2 

158.6 

120.2 

155.9 

166.5 

30% 

106 

105.6 

111.2 

107.1 

108.5 

124.1 

108.5 

96 

85.2 

75.6 

94.3 

94.9 

70.9 

108.8 

121.4 

133.9 

101.5 

129.7 

137.9 

40% 

60.4 

60.2 

65.3 

59.9 

59.5 

78.8 

60.2 

59.2 

55.4 

46.3 

59.9 

59.1 

45.4 

60.3 

72.3 

80 

59 

82.9 

88.2 

50% 

40 

40 

44.8 

40 

39.3 

56.1 

40.8 

40.2 

39.5 

30.2 

40.5 

42.3 

29.7 

40.2 

50.4 

53.7 

39.4 

55.9 

59.9 

Case  2 

60W0 

60C5 

60C10 

60C5H5 

60C7H7 

60C5T5 

60C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

537.7 

533.3 

539.4 

534.5 

532.3 

537.3 

536.5 

535.6 

535.4 

544.8 

538.6 

546.1 

552.6 

530 

523.7 

516.8 

505.1 

519.4 

508.1 

20% 

209.1 

207.1 

207.3 

207.3 

203.7 

218.2 

206 

197.9 

194.4 

190 

197.4 

192.6 

185.4 

199.1 

216.4 

224.7 

205 

216.6 

236.9 

25% 

149.2 

149.1 

154.8 

150.4 

148.9 

167.1 

154.6 

147.1 

135.9 

122.1 

143.9 

133.2 

115.5 

155.8 

173.7 

187.1 

155.2 

184 

205.4 

30% 

114.3 

114 

119.6 

114.2 

114.5 

130.9 

116.3 

111.7 

104.3 

89.7 

111.2 

104.4 

84.1 

119.3 

138.2 

147.4 

117.4 

145.3 

160 

40% 

65.6 

65.4 

70.9 

65.5 

65.3 

83.2 

67 

67 

63.3 

52.3 

66.1 

65.3 

47.8 

66.4 

80.4 

85.3 

65.9 

88.2 

95.3 

50% 

47.6 

47.4 

52.6 

47.2 

47.3 

64.9 

47.9 

47.1 

46.3 

35.9 

45.4 

49.8 

32.6 

49.5 

59.3 

62.4 

49 

65.1 

69.2 

Case  3 

60W0 

60C5 

60C10 

60C5H5 

60C7H7 

60C5T5 

60C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

373.2 

372.5 

376.2 

374.3 

369.4 

383.5 

391.7 

379.1 

347.4 

378.8 

380 

336.3 

381.7 

363.7 

385.5 

384 

361.1 

393.4 

392.6 

20% 

337.7 

339.8 

337.4 

349.8 

353.6 

345 

372.7 

368 

369.3 

398.5 

370.1 

384.1 

404.3 

332.9 

350 

336.8 

342.3 

353 

346 

25% 

175.7 

176.2 

179.7 

178.3 

176.6 

193.2 

182.2 

180.8 

170.4 

168.3 

174.3 

168.7 

167.6 

182.6 

207.4 

203.8 

181.4 

219.2 

213.2 

30% 

140.2 

138.6 

141.2 

138.8 

140.2 

162 

149.4 

143.4 

138.3 

135.3 

142.6 

138.8 

136.3 

155.4 

169.3 

170.9 

144.2 

181.8 

162.3 

40% 

110.7 

111.4 

115.8 

115.2 

116.2 

126 

113.8 

103.4 

98.6 

98.2 

104.6 

104.4 

97.6 

115 

122.8 

127 

113.6 

122.1 

128.9 

50% 

54.9 

55.4 

59.2 

57 

55.8 

69.5 

58.7 

63 

65.6 

61.1 

62.2 

68.6 

60.8 

57.5 

63.7 

62.6 

58.4 

68 

66.3 

Case  4 

60W0 

60C5 

60C10 

60C5H5 

60C7H7 

60C5T5 

60C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

242.9 

242.8 

248.2 

245 

239.8 

262 

260.3 

239 

233.6 

244.7 

241.7 

235.1 

256.7 

249.5 

259.6 

263.6 

247.6 

272.5 

265.4 

20% 

137.4 

139.4 

143.1 

140.8 

138.3 

159.6 

147.7 

137 

132.8 

116.5 

136.6 

133.3 

116.8 

140.6 

166.4 

178.3 

140.6 

173.6 

189.3 

25% 

114.3 

114 

122.3 

112 

109.9 

137.6 

118.3 

113.2 

115.6 

89.3 

112.4 

118.1 

82.1 

107.6 

132.8 

148.6 

109.1 

141.1 

161.8 

30% 

88.7 

89 

96.7 

90.1 

88.3 

112.2 

96.2 

89.4 

88.1 

71.1 

89.1 

91.9 

70.2 

86.5 

102.7 

118.3 

84.5 

109.9 

127.4 

40% 

53.7 

53.7 

61.4 

53.7 

52.6 

75 

57.2 

53.3 

52.1 

39.2 

53 

56.2 

38.9 

51.5 

66.3 

73.5 

52.3 

73.2 

83 

50% 

38.8 

38.9 

42.5 

38.6 

38.2 

52 

40.7 

36.7 

36.4 

26.6 

36.6 

39 

25.6 

37.8 

46.8 

52.2 

37.8 

52 

58.3 

From  Tables  V-VII,  an  interesting  observation  is  that  for  the  situations  with  headwind  and  very  short  observation  time, 
e.g.,  C5H5  and  C7H7  with  small  OTs  (e.g.,  OT=15%  or  20%),  the  g.w.i.  IPP  RMS  performance  is  counter-intuitively  worse 
than  the  n.w.i.  IPP  RMS  performance.  This  is  because  of  the  imperfect  selection  of  which  is  the  best  IMM  in  the  MIMM 


estimator  when  it  is  based  on  a  very  short  observation  time.  More  clearly,  as  shown  in  Fig.  5  (the  trajectory  60C7H7  in  Case 
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TABLE  VI:  IPP  RMS  errors  (in  m,  100  MC  runs)  for  various  OTs  (shown  in  the  first  column)  and  various  sensor  accuracies, 

caliber  81mm:  italic,  see  footnote  8;  bold,  see  footnote  9. 


Case  1 

81W0 

81C5 

81C10 

81C5H5 

81C7H7 

81C5T5 

81C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

420.8 

420.1 

424.1 

422.2 

420 

438.7 

427.4 

401.8 

423.3 

428.6 

389.3 

436.2 

399.1 

421.1 

382.8 

414.3 

420.3 

384.9 

415 

20% 

432.4 

428.4 

438.6 

429 

428 

445.9 

430.9 

423.4 

401.8 

384.1 

418.6 

393.2 

361.2 

428.7 

456.5 

481.9 

430.1 

475.5 

497.6 

25% 

288.9 

294.4 

294.5 

292.3 

292.7 

300.6 

296.3 

292.3 

260.8 

254.7 

291.1 

257.4 

240.3 

282.7 

323.6 

333 

285.2 

343.8 

352.6 

30% 

188 

186 

193.3 

186.8 

185 

204.1 

190 

190.2 

181.2 

158 

196.6 

183.3 

148.5 

189.3 

206 

230.1 

191.5 

221.3 

246 

40% 

101.1 

100.5 

104.8 

101 

100.8 

118.5 

103.4 

99.4 

96.7 

80.4 

99 

96.8 

72.7 

99.3 

114.2 

127.1 

101 

122.4 

141.6 

50% 

48.7 

48.8 

51.8 

48.6 

49.7 

62.2 

51.6 

48.5 

47.9 

40.1 

47.7 

49.7 

36.8 

48.5 

58.6 

63.7 

49.5 

64.5 

71.9 

Case  2 

81W0 

81C5 

81C10 

81C5H5 

81C7H7 

81C5T5 

81C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

497.5 

492 

492.2 

491.4 

482.9 

501.6 

484.2 

490.5 

520.3 

500.4 

511.8 

504.4 

528.4 

491.9 

526.6 

494.1 

486.7 

534.8 

500.7 

20% 

450.5 

449.6 

451.7 

448.9 

448.9 

459.7 

448.8 

451.7 

416.9 

398.2 

441.6 

401.2 

383 

444.3 

467.6 

486.7 

449.9 

495.5 

510.2 

25% 

315.7 

313.6 

320 

317 

313.1 

328.4 

318.4 

318.8 

295.6 

283 

317.8 

285.9 

269.3 

317.1 

354.3 

359.6 

320 

376.9 

377.2 

30% 

209.5 

209 

217.1 

209.8 

208.4 

230.2 

211.6 

209.3 

204.6 

175.3 

217.2 

214.8 

174.7 

212.5 

239.6 

251.8 

213.7 

257.8 

267.3 

40% 

133.8 

135.2 

140 

135.4 

134.9 

151.3 

137.1 

135.3 

130 

111.5 

134.9 

130.9 

103.3 

131 

150.9 

157.7 

129.3 

159.9 

166.4 

50% 

68.9 

69.1 

73 

68.6 

72.7 

87.5 

73.6 

72.1 

69.7 

57.9 

71.8 

71.7 

53.1 

69.1 

80.1 

86.1 

69.6 

87.7 

94.7 

Case  3 

81W0 

81C5 

81C10 

81C5H5 

81C7H7 

81C5T5 

81C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

665.2 

665 

667.1 

665.9 

652.6 

675.7 

6J9.2 

664.8 

661.7 

643.9 

690.8 

703.7 

668.2 

646.5 

673.5 

667.4 

651.7 

687.6 

667.4 

20% 

479.9 

479.6 

480.9 

479.7 

485 

490.7 

490.9 

478.3 

473.4 

456 

472.8 

471.7 

445.3 

476.8 

473.2 

504.8 

476.6 

501.9 

516.7 

25% 

344.9 

345.1 

349.3 

346.4 

348.2 

362.3 

352.3 

351.4 

327.9 

331 

348.3 

320.9 

325.5 

347.4 

375.1 

365.2 

346.9 

372.1 

376 

30% 

213.5 

213.7 

216.9 

212 

217.5 

229.5 

218.2 

212.7 

201.5 

190.4 

223.4 

203.3 

208.9 

217.9 

240.9 

238.4 

219.1 

248.8 

252.2 

40% 

131.2 

127.2 

136.6 

126.9 

126.9 

144.4 

131.7 

131.4 

130.7 

124.4 

128.4 

134.1 

129.7 

131.6 

145.5 

139.4 

129.9 

143 

145.2 

50% 

66.7 

67 

70.4 

67 

63.9 

74.8 

67.1 

62.5 

63 

60.4 

64 

65.3 

60 

66.8 

71.1 

76.4 

67.7 

77.8 

83 

Case  4 

81W0 

81C5 

81C10 

81C5H5 

81C7H7 

81C5T5 

81C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

465.4 

465 

468.2 

465.2 

464.9 

479.5 

467.3 

464.9 

444.7 

395.7 

445.1 

451.4 

375.8 

464.7 

469.2 

520.1 

463 

483 

549.7 

20% 

322.4 

322 

332.9 

322.3 

321.5 

342.7 

323.6 

324.1 

294.5 

268 

323.2 

286.4 

256.3 

317.9 

362.1 

368.7 

301.1 

367.5 

375.8 

25% 

199.8 

202.1 

208.2 

203.1 

202 

223.4 

206.9 

200.4 

176.8 

165 

198.5 

174.7 

154 

201.8 

235 

242.5 

205.2 

252.5 

261.4 

30% 

123.9 

122.4 

129.1 

123.2 

121.3 

145.6 

130.2 

121.1 

120.6 

95.8 

124.3 

138 

94.3 

122.6 

140.2 

166.6 

121.2 

151.4 

184.8 

40% 

67.7 

67.7 

73.5 

68.8 

67.7 

85.2 

72.6 

67.4 

65.6 

52.7 

67.5 

68.2 

50.9 

67.8 

80.8 

91.8 

65.7 

87.7 

101.6 

50% 

33.5 

33.2 

37.3 

34.7 

36.2 

48.6 

40.8 

36.2 

35.5 

28.1 

35.9 

37.6 

28.3 

33.8 

42.7 

49.2 

34 

48.2 

58.1 

3  and  OT=15%),  because  of  this  imperfect  decision  of  the  best  IMM,  each  MC  IPP  “cloud  of  points”  (for  g.w.i.  and  n.w.i.) 
is  approximately  separated  into  3  different  clusters,  designated  as  the  “cluster  A”,  “cluster  B”  and  “cluster  C”.  It  can  be  seen, 
compared  with  the  g.w.i.  clusters,  that  the  n.w.i.  clusters  have  a  clear  drift  caused  by  the  wind  effect.  By  coincidence,  this 


drift  happens  to  drive  “cluster  A  (n.w.i.)”  that  contains  the  most  points  of  the  MC  IPP  closer  to  the  true  impact  point.  This 
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TABLE  VII:  IPP  RMS  errors  (in  m,  100  MC  runs)  for  various  OTs  (shown  in  the  first  column)  and  various  sensor 
accuracies,  caliber  120  mm:  italic,  see  footnote  8;  bold,  see  footnote  9. 


Case  1 

120W0 

120C5 

120C10 

120C5H5 

120C7H7 

120C5T5 

120C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

718.2 

718.2 

720 

718.1 

717.8 

714.7 

718.8 

749.5 

653.1 

714.2 

753.2 

652.5 

698.9 

684.6 

760.6 

723.1 

681.3 

779.5 

745.6 

20% 

547.4 

546.7 

548.5 

547.4 

539.6 

555.6 

544.5 

576.7 

562.3 

534.9 

572.4 

575.6 

477.1 

581.2 

571.5 

651.9 

591.7 

585.5 

699.5 

25% 

488.9 

492.6 

492 

492.9 

491.9 

497.4 

497.9 

488.5 

474 

416.3 

484.2 

465.9 

391.2 

502.9 

521.7 

552.4 

494.8 

513.7 

578 

30% 

308.6 

308.4 

315.5 

308.4 

307.8 

330.4 

308.9 

306.4 

292.7 

261.8 

302 

284.7 

243.5 

313 

334 

354.3 

316 

337.8 

368.4 

40% 

117.5 

115 

122.6 

115.8 

114.8 

132.6 

117 

114.6 

116.6 

100.4 

115.4 

118 

94.9 

121.5 

130.1 

142.2 

118.2 

128.7 

150.7 

50% 

55.1 

54.9 

58.8 

56.6 

56 

67.2 

59.8 

52.4 

55.4 

52.9 

51.4 

57.5 

54 

56 

60.6 

62.6 

56.7 

64.1 

66.8 

Case  2 

120W0 

120C5 

120C10 

120C5H5 

120C7H7 

120C5T5 

120C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

779.2 

776.3 

790.9 

775.8 

782.2 

789 

777.1 

769.3 

748.9 

766.6 

786 

725.8 

756.8 

761.2 

797.5 

773.3 

744.9 

788.3 

795.7 

20% 

613.3 

613.3 

618 

612.3 

616 

622.6 

614.8 

597.5 

593.1 

535.7 

607.2 

571.8 

509.9 

624.1 

631.1 

693.5 

621.4 

653.9 

725 

25% 

489.8 

494.1 

492.5 

493.5 

493.6 

497.9 

492.9 

493.2 

480.9 

450 

497.5 

473.7 

419.3 

493.9 

500.7 

557.1 

475.7 

486.2 

551.6 

30% 

346.7 

348.2 

349.8 

347.7 

347.9 

354.3 

345.2 

362 

346.4 

317.5 

359.7 

344.5 

299.6 

345.9 

368.2 

387.8 

342 

381.2 

414.3 

40% 

144.5 

142.9 

148.8 

143.5 

142.3 

158.3 

144.4 

144.9 

148.7 

126.7 

144.9 

149.1 

119.6 

152.3 

159.1 

179.4 

156.4 

174.1 

197.2 

50% 

63 

61.5 

66.7 

62 

61.1 

75.7 

63.4 

61.6 

63.4 

53.7 

59.8 

64.8 

50.4 

64.7 

69.1 

79.6 

63.1 

73.8 

84.9 

Case  3 

120W0 

120C5 

120C10 

120C5H5 

120C7H7 

120C5T5 

120C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

887.6 

879.9 

885.9 

879.6 

879.2 

889 

879.3 

874.7 

859.6 

833.7 

870.8 

867.2 

817.1 

883.2 

894.5 

918.5 

878.7 

901.2 

936.4 

20% 

650.1 

651.2 

654.4 

652.9 

648.6 

659.1 

643.8 

634.1 

611.4 

612.4 

642.2 

608.4 

607.2 

639.8 

686.1 

673.8 

642.7 

697.9 

685.1 

25% 

446.3 

446.4 

443.7 

433.6 

431.7 

445 

431.3 

449.1 

397.7 

433.1 

437.6 

395.9 

425.4 

426.6 

465.4 

450.9 

398.9 

422.2 

452.6 

30% 

279.3 

271.4 

281.5 

272.5 

270.1 

302.4 

275.9 

284.4 

277.7 

273.1 

285.7 

275.8 

265.9 

272.5 

284.1 

294.4 

272.4 

286.6 

314.3 

40% 

130.5 

130.8 

135.6 

131 

127.1 

146.5 

133 

136.7 

135.4 

132.3 

134.3 

138.6 

131.7 

129.1 

139.5 

132.4 

133.3 

151.4 

144 

50% 

64.9 

65.1 

68.1 

67.7 

64.4 

75 

71 

64.3 

67.3 

65.5 

64.1 

70.2 

66.6 

66.7 

70.3 

68.1 

62.6 

68.1 

71.5 

Case  4 

120W0 

120C5 

120C10 

120C5H5 

120C7H7 

120C5T5 

120C7T7 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

g.w.i. 

n.w.i. 

e.w.i. 

15% 

462.7 

462.4 

462.2 

477.3 

476.6 

468.3 

482.2 

480.2 

445.6 

426.4 

464.8 

456.8 

453.1 

464.1 

534.3 

533.8 

474.7 

528.1 

590.9 

20% 

433.9 

433.3 

433.8 

434.1 

430.2 

448.5 

434.6 

413.3 

420.6 

350.8 

413.8 

445.6 

306.8 

444 

430.8 

506 

456 

437.5 

543.7 

25% 

337.4 

336.6 

340.7 

336.9 

335.5 

353.1 

339.8 

323.9 

318.7 

276 

321.8 

321.2 

260 

350.1 

360.3 

394.8 

329.5 

368.2 

398.2 

30% 

198.1 

198 

206.4 

196.7 

196 

213.1 

200.7 

177.9 

170.4 

151.8 

178.9 

172.2 

143.1 

207.9 

221.4 

246.6 

211 

235.1 

268.4 

40% 

62.2 

62 

65.1 

62.9 

61.4 

76 

67.6 

61.8 

65.8 

59.7 

59.9 

69.2 

64.2 

64.5 

68.9 

80.2 

67.7 

77.8 

93 

50% 

29.5 

29.4 

32.4 

30.9 

29.4 

40.7 

35 

30.5 

33 

32.5 

29.8 

36.2 

34.9 

30.2 

34.6 

38.6 

29.7 

37.9 

43.2 

results  in  the  g.w.i.  IPP  RMS  errors  being  worse  than  the  corresponding  n.w.i.  IPP  RMS  errors  in  some  cases  in  Tables  V-VII. 
However,  as  the  observation  time  increases,  the  MIMM  estimator  can  gradually  make  the  correct  decision  to  find  the  best 
IMM  estimator  (the  one  with  the  most  suitable  initial  drag  coefficient  estimate)  and  then,  as  we  expect,  the  g.w.i.  IPP  results 


are  better  than  the  corresponding  n.w.i.  IPP  results.  Given  a  very  short  OT  in  a  strong  wind  environment  (e.w.i.  vs.  g.w.i.),  the 
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bias  in  the  erroneous  wind  information  sometimes  exacerbates  and  sometimes  improves  the  IPP  results.  Generally,  the  e.w.i. 
IPP  performance  is  between  the  g.w.i.  and  n.w.i.  cases  when  the  MIMM-selected  drag  coefficient  estimate  matches  the  truth. 

3400 

3200 

3000 

2800 
c 
o 

w 

?-  2600 

2400 

2200 

2000 

Fig.  5:  The  clusters  of  MC  IPP  “clould  of  points”  caused  by  imperfect  decision  of  the  best  IMM  in  MIMM  estimator. 
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B.  IPP  Uncertainty  Ellipse 

With  various  OTs  and  sensor  accuracies,  the  IPP  uncertainty  ellipses  using  (22),  centered  at  the  true  impact  point,  with  the 
corresponding  MC  IPP  “cloud  of  points”  (for  both  g.w.i.  and  n.w.i.)  for  various  trajectories  in  the  presence  of  wind  are  shown 
in  Figs.  6-8. 

From  Fig.  6,  it  can  be  seen  that  as  the  OT  increases,  more  IPP  points  fall  in  the  corresponding  uncertainty  ellipses  and  the 
sizes  of  the  uncertainly  ellipses  gradually  shrink.  The  total  observation  time  is  crucial  for  the  evaluation:  given  the  OT  is  large 
enough  (50%),  most  of  the  MC  IPP  points  fall  within  the  uncertainty  ellipses  with  acceptable  IPP  RMS  errors.  The  outliers 
are  due  to  the  forced  early  decision  when  the  OT  is  small. 

From  Fig.  7,  it  can  be  seen  that  as  the  sensor  accuracy  worsens,  the  corresponding  IPP  uncertainty  ellipses  become  wider 
because  of  the  poor  angle  accuracy  or  become  longer  because  of  the  poor  range  accuracy,  and  fewer  MC  IPP  points  fall  into 
them.  Particularly,  the  larger  sensor  errors  in  angles  (which  cause  more  uncertainty  along  the  minor  axis  of  the  uncertainty 
ellipse  in  the  case  considered)  result  in  the  MC  IPP  points  being  more  scattered  along  the  minor  axis  of  the  corresponding 


IPP  uncertainty  ellipse;  the  larger  sensor  errors  in  range  (which  cause  more  uncertainty  along  the  major  axis  of  the  uncertainty 
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ellipse  in  the  case  considered)  result  in  the  MC  IPP  points  being  more  scattered  along  the  major  axis  of  the  corresponding  IPP 
uncertainty  ellipse.  Overall,  the  sensor  errors  along  the  projectile  travel  direction  have  more  impact  on  the  IPP  performance. 

From  Fig.  8,  comparing  the  n.w.i.  IPP  clouds  with  smaller  wind  (“C5’7“C5H5”/“C5T5”)  to  those  n.w.i.  IPP  clouds  with 
stronger  wind  (“C10’7“C7H7”/“C7T7”),  respectively,  it  is  obvious  that  the  stronger  the  wind  is,  the  more  significant  the 
deviations  of  the  n.w.i.  IPP  clouds  are.  When  both  the  cross  wind  and  the  range  wind  were  present,  the  wind-induced  drift 
of  the  n.w.i.  IPP  clouds  is  compounded;  without  accounting  for  the  cross  wind,  the  n.w.i.  IPP  cloud  deviates  to  one  side 
(the  right  side  in  the  cases  discussed)  of  the  major  axis  of  the  uncertainty  ellipses  (obtained  by  properly  accounting  for  the 
wind  effect  in  the  dynamic  model)  while  the  lack  of  accounting  for  the  range  wind  causes  the  n.w.i.  IPP  cloud  to  be  “pushed 
backward/forward”  causing  range  under-Zover-  prediction. 

C.  Consistency  Evaluation 

The  consistency  of  the  MIMM  estimator  is  very  important  in  the  IPP  application.  It  is  always  desired  that  the  predicted 
impact  point  falls  into  the  uncertainty  ellipse  centered  at  the  true  impact  point  (or  the  equivalent  converse).  This  can  be  achieved 
by  increasing  the  process  noise  intensities,  which  results  in  an  ellipse  with  its  major-  and  minor-  axes  large  enough  to  include 
the  predicted  impact  point.  However,  this  solution  is  undesirable.  A  consistency  test  will  help  to  find  the  process  noise  setting 
that  gives  the  uncertainty  ellipse  compatible  with  the  actual  errors  [1]. 

The  consistency  of  the  MIMM  estimator  is  examined  using  both  the  normalized  estimation  error  squared  (NEES),  which 
is  preferable  for  Monte  Carlo  runs  when  the  truth  is  available  (off-line  simulations),  and  the  normalized  innovation  squared 
(NIS).  The  latter  is  the  only  one  that  can  be  used  in  real  time  testing  [1]. 

Eig.  9  shows,  for  the  best  IMM  estimator  selected  from  the  MIMM  estimator  with  OT=50%  and  various  sensor  accuracies, 
the  NEES  consistency  in  position  and  velocity  and  the  NIS.  Note  that  the  NIS  is  evaluated  for  TM  before  the  estimated 
BoP  and  for  BM  afterward  [15].  It  can  be  seen  that  as  the  sensor  errors  increase,  the  IMM  estimator  selected  becomes  less 
consistent.  In  particular,  the  larger  the  sensor’s  range  error  (which  causes  more  uncertainty  along  the  travel  direction  of  the 
projectiles  in  present  cases),  the  more  significant  is  the  inconsistency. 

VI.  Summary  and  Conclusions 

The  MIMM  estimator  developed  in  [15]  has  been  extended  to  account  for  the  wind  effect.  The  IPP  performance  for  various 
total  observation  times  and  sensor  accuracies  have  been  investigated  in  detail.  The  wind  effect  can  be  fully  compensated  if 
the  wind  information  is  available;  it  also  can  be  mitigated  by  increasing  the  total  observation  time  if  no  wind  information  is 
provided.  The  total  observation  time  is  key  for  IPP  performance:  for  example,  given  the  observations  all  the  way  up  to  the 
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apogee  one  can  always  expect  a  very  good  IPP  performance,  but  this  is  clearly  not  desirable  if  countermeasures  are  to  be  taken 
against  the  projectile.  The  -point  initialization,  designed  for  initializing  the  key  thrust  component  as  well  as  the  kinematic 
components,  is  considered  a  good  method  to  quicken  the  correct  mode  convergence  in  each  IMM  estimator  and  alleviate  the 
estimation  ambiguity  between  the  drag  coefficient  and  the  thrust  components.  In  the  cases  studied  here,  an  OT  of  20%  or 
more  and  sensor  accuracies  under  10  m  in  range  and  5mrad  in  angles  are  required  to  overcome  the  marginal  observability 
(ambiguity)  of  the  drag  and  thrust  and  then  achieve  an  acceptable  IPP  performance,  with  the  IPP  RMS  error  no  more  than 
10%  of  the  projectile  range.  The  sensor  errors  that  cause  more  uncertainty  along  the  travel  direction  of  the  projectiles  have 
more  impact  on  the  IPP  performance.  In  the  presence  of  wind,  the  short  observation  time  and  the  limited  sensor  accuracy  are 
more  critical  since  they  could  cause  confusion  in  selecting  the  “correct”  best  IMM  estimator  (which  should  be  the  one  with 
the  most  suitable  initial  drag  coefficient  estimate,  closest  to  truth)  from  the  MIMM  estimator,  and  then  lead  to  a  degraded  IPP 
performance. 


Appendix  A 

Cubic  Spline  Interpolation 

As  shown  in  Fig.  10,  the  cubic  spline  curve  of  the  Mach  number-dependent  drag  coefficient  multiplier  is  obtained  by 
interpolating  the  selected  representative  points  shown  in  Table  VIII. 

TABLE  VIII:  Selected  representative  points  for  cubic  spline  interpolation 


Speed  (m/s) 

1 

50 

too 

150 

190 

230 

280 

285 

295 

310 

350 

Normalized  drag  coefficient 

1 

I.OI 

1.015 

1.02 

1.03 

1.035 

1.1 

1.2 

1.3 

1.4 

1.59 

Speed  (m/s) 

400 

440 

500 

550 

590 

650 

750 

800 

850 

1000 

1000  + 

Normalized  drag  coefficient 

1.69 

1.71 

1.68 

1.55 

1.4 

1.35 

1.33 

1.32 

1.31 

1.3 

1.3 
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Fig.  6;  IPP  uncertainty  ellipse  and  MC  IPP  cloud  for  various  OTs,  trajectory  81C5T5,  Case  1,  100  MC  runs 
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Fig.  7:  IPP  uncertainty  ellipse  and  MC  IPP  cloud  for  various  sensor  accuracies,  trajectory  60C10,  OT=25%,  100  MC  runs 
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Fig.  8;  The  IPP  uncertainty  ellipse  and  MC  IPP  cloud  for  various  120mm  trajectories,  OT=30%  100  MC  runs 
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Fig.  9;  Consistency  test  (NEES  in  pos.  and  vel.  and  NIS),  trajectory  81C7T7,  OT=50%,  100-MC  run 
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Fig.  10:  The  cubic  spline  approximation  (for  a  “sharp  nose”  projectile  [17]). 


